# -*- coding: utf-8 -*-
"""
Created on Fri Jan 01 15:38:55 2016

@author: HP
"""
import math

def Gauss(longitude,latitude):
    ProjNo = 0
    ZoneWide = 6
    iPI = 0.0174532925199433  #3.1415926535898/180.0
    a=float(6378140)
    b=float(6356755)
    ProjNo = (int)(longitude / ZoneWide) 
    longitude0 = ProjNo * ZoneWide + ZoneWide / 2
    longitude0 = longitude0 * iPI
    #latitude0 = 0
    longitude1 = longitude * iPI  #经度转换为弧度
    latitude1 = latitude * iPI   #纬度转换为弧度
    e2=(a**2-b**2)/(a**2)   #第一第二偏心率的平方
    ee=float((a**2-b**2)/(b**2))

    NN=a/math.sqrt(1.0-e2*math.sin(latitude1)*math.sin(latitude1))
    T=math.tan(latitude1)*math.tan(latitude1)
    C=ee*math.cos(latitude1)*math.cos(latitude1)
    A=(longitude1-longitude0)*math.cos(latitude1)
    M=a*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude1-(3*e2/8+3*e2*e2/32+45*e2*e2*e2/1024)*math.sin(2*latitude1)+(15*e2*e2/256+45*e2*e2*e2/1024)*math.sin(4*latitude1)-(35*e2*e2*e2/3072)*math.sin(6*latitude1))
    xval = NN*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120)
    yval = M+NN*math.tan(latitude1)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24+(61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720)
    X0 = 1000000L*(ProjNo+1)+500000L
    Y0 = 0
    xval = xval+X0
    yval = yval+Y0
    return xval,yval
    
    
print Gauss(114,50)